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■^^ ' We give a well posed initial value fonnulation of the Baumgarte-Shapiro-Shibata-Nakamura form 

f^ , of Einstein's equations with gauge conditions given by a Bona-Masso like slicing condition for the 

f^ ■ lapse and a frozen shift. This is achieved by introducing extra variables and recasting the evolution 

^S| ' equations into a first order symmetric hyperbolic system. We also consider the presence of artificial 

boundaries and derive a set of boundary conditions that guarantee that the resulting initial-boundary 
value problem is well posed, though not necessarily compatible with the constraints. In the case of 
dynamical gauge conditions for the lapse and shift we obtain a class of evolution equations which 
are strongly hyperbolic and so yield well posed initial value formulations. 

^ . I. MOTIVATION 

<^ . 

^D . Most numerical evolutions of Einstein's field equations try to approximate solutions on a generically infinite 

^^ (non-compact) 3-space by computations on a truncated finite (compact) domain. For this, artificial boundaries and 

corresponding boundary conditions have to be introduced. Mathematically, this immediately poses the question of 
well-posedness of the initial-boundary-value problem (IBVP) for the evolution equations and compatibility with the 

(f^ ] constraints. In addition, governed by causality, the solution on the finite domain is very likely to differ from the 
solution on the infinite domain, after disturbances from the boundaries enter the computational domain. This makes 
the choice of the boundary condition crucial for the physical interpretation of the results, especially if one thinks of 
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I . integrated quantities like masses, charges and momenta. Physically, it could be even argued that ultimately only 









estimates of the deviation of the numerical solution from the actual solution are significant and that the mathematical 
concepts of well-posedness of the IBVP including compatibility with the constraints and avoidance of reflections from 
artificial boundaries are only steps towards achieving this goal. 



5_j ■ Removing the influence of the boundaries could be achieved by enlarging the computational domain to a size 
Ci ' such that, according to causality, disturbances from boundaries cannot have reached the domain of physical interest. 
Note that this requires knowledge of the causality structure of the spacetime, which presupposes estimates on the 
solution to be found. In addition, increasing the size of the computational domain goes at the cost of resolution, 
because of finite computational resources, which is particularly restricting in 3— dimensional problems although 
this problem can be alleviated by using adaptive or fixed mesh refinement techniques [J]. So, both restrictions in 
computational resources and the demand for higher resolution, lead us to attempt to minimize the influence of 
artificial boundaries on the numerical solution. This is tried by so called 'outgoing boundary conditions' meant to 
make those boundaries to appear as 'transparent' as possible. For instance, one such approach is given by Endquist 
and Majda |3| using a hierarchy of conditions which gradually decrease reflections at the boundary. Alternative 
approaches are the methods of characteristic ^J or perturbative matching 4] . See [3, |6| for an approach trying to 
avoid the introduction of artificial boundaries altogether by a suitable compactification of spacetime. 

II. INTRODUCTION 

In this article we analyze the IBVP of the Baumgarte-Shapiro-Shibata-Nakamura (BSSN) [3, |a| formulation of 
Einstein's vacuum equations, which is currently used by several groups in numerical relativity with applications to 
the binary black hole and binary neutron star problem, see |9] for a review. Since the BSSN equations are first 
order in time, but mixed first/second order in space, their type (elliptic, parabolic, hyperbolic or mixed) is a priori 



not clear. Here, we analyze the well posedness of their (nonlinear) Cauchy problem with and without boundaries. 
Well posedness means that the Cauchy problem has a unique solution local in time and that the solution depends 
continuously on the initial data. The last property is important in view of obtaining convergent discretizations 
since in general numerical simulations introduce small errors in the initial data. If violated, this can lead to er- 
rors at a later time which grow exponentially with increasing resolution (see llOllllj for examples of this phenomenon). 

We find that the BSSN system with a large family of gauge conditions for the lapse, including Bona-Masso 
like slicing conditions, and an a priori specified shift yields a well posed initial value problem. This is achieved by 
introducing extra fields that make it possible to recast the system into a first order quasi-linear symmetric hyperbolic 
form for which standard well posedness results are known [IJ. The introduction of extra fields brings additional 
constraints, and the original BSSN system and the first order symmetric hyperbolic system derived in this article are 
only equivalent if these constraints are satisfied. However, we show that the associated constraint variables obey a 
closed evolution system that is independent of the other constraints. This means that the additional constraints are 
satisfied everywhere at later times if satisfied initially, even if the other constraints are violated. This implies that 
the (original) BSSN system is well posed; in particular, unique solutions local in time exist, and depend continuously 
on the initial data. 

Our first order symmetric hyperbolic reduction also facilitates the analysis of characteristic modes which is 
particularly useful when constructing boundary conditions. Here we construct maximally dissipative boundary con- 
ditions that guarantee the well posedness of the resulting IBVP [13 ■ Although in general these boundary conditions 
are not compatible with the constraints, they are consistent with the evolution equations and constitute a first step 
towards improving numerical evolutions of the BSSN system. In particular, the present analysis offers the possibil- 
ity to construct constraint- preserving boundary conditions 14] in the linearized case, following the lines of |l5lll6lll7| . 

The techniques used in this article are the same used in [13 where well posedness of the BSSN system with 
an explicitly given shift and an algebraic gauge condition is found by considering an auxiliary first order system. A 
different technique which makes use of pseudo-differential calculus has recently been applied in order to show well 
posedness for a closely related formulation [13 . More recently, in [13, 113 a definition of symmetric hyperbolicity 
based on energy estimates for second order systems was presented which was verified for the BSSN system and the 
formulation in 19] for the case of an algebraic lapse and an explicitly given shift. Nevertheless the connection of 
their definition and existence of solutions is open. 

The remainder of this work is organized as follows. In section IIIII we review the BSSN equations, specify the 
gauge conditions we are considering and discuss the evolution system for the constraint variables. In section IIVI 
we introduce extra fields and derive a first order symmetric hyperbolic system that refiects the dynamics of the 
original BSSN system. The characteristic fields with nontrivial speeds are constructed in section and are used 
to write down maximally dissipative boundary conditions. In section I VII we find using pseudo-differential calculus 
that the BSSN system with a "iiT-driver" and a "Gamma-freezing" condition as defined in |21| but with a different 
time coordinate is strongly hyperbolic according to the definition in |l9l |22| | and so yields a well posed initial value 
formulation. Conclusions are drawn in section rVIII 



III. THE BSSN EQUATIONS 

Since our results depend crucially on the principal part of the equations, we write down the BSSN system explicitly 
in this section. The system of equations is the one that has been used in Ref. [21| for numerical simulations, but it 
might differ from the one used by other groups. Decomposing the three metric and the extrinsic curvature according 
to 

7., = e''^, , (1) 

K^, = e^* (iy + I^^.k) , (2) 



where "/ij has unit determinant and K = ^^^ Kij is the mean curvature, the evolution equations are obtained from 

6^+6 



9o0 = -'^K+lduP^ (3) 
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do% = -2aAy -f 27fc(,9j)/?'' - -7y9fe/3^ (4) 
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+ aKA^j ~ 2aA,kA) + 2Ak{,dj)(3'' 



'^AjdkP'' 



— aS, 



-40 A 



dor = j^'dkdiP' + -f'd.dkP'^ + dki""' ■ d,p' - -dki''' ■ d,p' 
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2m ~ ^ _ _ ' 

\)dkA''' - —D'K + m{T\iA''^ + dA'^djCJ)) 



S\ 



(5) 



(6) 



(7) 



where we have mtroduced the operator do — dt ~ P-'dj. Here, all quantities with a tilde refer to the conformal three 
metric jij, and the latter is used in order to raise and lower their indices. In particular, Di and F ^ refer to the 
covariant derivative and the Christoffel symbols, respectively, with respect to jij. The expression [...] denotes the 
trace-less part (with respect to the metric 7^ ) of the expression inside the parentheses, and 



Rij - 



Rl 



= -\f'dkdi% + 7fe(.9,)f'= - f(y),9,7^'= + 7'^ (2f'5(,f,)fc, + fifH,) , 



-2D,Dj(f> - 2%jb''bk(j) + Ab,(f>bj(t) - A^,jb^cj)bk 



The parameter m, which was introduced in 23], controls how the momentum constraint is added to the evolution 
equations for the variable r\ The system in Ref. 21j corresponds to the choice m = 1. However, in order to obtain 
a first order symmetric hyperbolic reduction, we will see later that we need to to be a specific function of the lapse 
and the mean curvature. The source terms S, Sij and S"* are defined in terms of the four Ricci tensor, i?), , and the 
constraint variables 
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S' 



2amf^Mj~doCl 



(8) 

(9) 
(10) 

(11) 
(12) 
(13) 



The vacuum equations consist of the evolution equations H3I7|) with 5 = 0, 5*^ = 0, 5* = and the constraints H — 0, 



Mi = and Ci- 



0. 



Using the Bianchi identities, 2V^R\j,J, — V vR'^^^ = and imposing the evolution equations, it can be shown that 



the constraint variables obey the following propagation system: 

1 



doH = --D^a^Mj) 
a 

v3 



2a 



ae-4^i'^7fe,9,C^- + ^ifiI, 






2afnj'''Mi . 



(14) 

(15) 
(16) 



By introducing the further constraint variable Z,'' — djCp which satisfies d[iZ .J' = one can reduce Eqs. (|14I15I16|) 
to a first order symmetric hyperbolic system provided that m > 1/4. In the absence of boundaries, this implies that 
the constraints are preserved, i.e. trivial initial data for the constraints variables lead to zero constraint variables 
at later times as well. If time-like boundaries are present, the constraints are only preserved if suitable boundary 
conditions are specified. Such constraint-preserving boundary conditions are discussed in |l4Lll5lll6lll7|: but for the 
(nonlinear) BSSN system it is not yet understood if they lead to a well posed IBVP. 



In order to evolve the system (|3I7|I we have to specify conditions on the lapse a and the shift /3^ The simplest 
possibility is to set a = 1 (or any other fixed function) and /?* — 0. However, this leads to a formulation that is not 
strongly hyperbolic (this will follow from the results in Sect. IVII if we set the function f defined below in Eq. H21|l to 
zero). This can be avoided by "densitizing" the lapse. More generally, we can require [2j that the lapse 

a = a(0,2:^), (17) 

is a smooth strictly positive function of the conformal factor (or the determinant of the three metric) and spacetime 
coordinates with the restriction that a = (12a)^^9a/90 is strictly positive. Taking a time derivative of this, assuming 
that da/dt = and using Eq. ^ we obtain 

^a = -2aV (k - -DhP""] , (18) 

at \ a J 



which is the modification of the Bona-Masso condition [25| proposed in |26l |27| . The advantage of this gauge is that 
it is compatible with a time-independent lapse in a time slicing that is adapted to stationarity if dt is a Killing field. 
It follows from the calculations of Ref. [Tg that in this case the BSSN system is strongly hyperbolic if one chooses 
ra > 1/4 and symmetric hyperbolic if the parameter ra is adjusted such that Am = 6cr + 1 with a > 1/2 JSQ- We 
mention here that the special case a = e^'^Q{x^), where Q{x^) is an a priori specified function, has been observed to 
lead to more stable numerical evolutions of single black hole with the BSSN system p^. 

Here, we are interested in live gauge conditions which allow lapse and shift to react on changes of the fields. Such 
conditions can be useful, for instance, to evade singularities. In this article, we consider two cases of gauge conditions: 

(a) The following evolution equation for the lapse 

doa^-aF{a,K,x^'), (19) 

where _F is a smooth function of a, K and x^ with the restriction that 

I dF 

This condition generalizes the Bona-Masso gauges. The shift is frozen, that is, assumed to be an a priori 
specified function of spacetime. Symmetric hyperbolic formulations of the vacuum field equations with these 
gauge conditions were obtained in |24j . 

(b) The gauge conditions of Ref. [21| which, for the lapse, require the "hyperbolic i^-driver" condition 

doa = -a^f{a,cl>,xn{K ~ Ko(xn), (21) 

where the function f{a, cj), x^) is smooth and strictly positive, and Kq(x^) is an arbitrary smooth function. For 
the shift, the "hyperbolic Gamma driver" [21| type condition 

do(3' = a^G{a,<j>,x'')B\ (22) 

doB' = e-'^'''H{a,<j),x'')dor-T]{B\a,x'') (23) 

is imposed, where G{a, (j), x^) and H(a^ 4>, x^) are smooth, strictly positive function, and ri{B^ , a, x^) is a smooth 
function. Notice that Eq. ()21|l is a special case of Eq. (|19|l . Note also that the conditions H21I22I23|) differ from 
the ones considered in J2j| by the replacement dt i— > i9o which simplifies the analysis in the present article. 

In the next section, we show that the gauge conditions (a) lead to a well posed initial value problem provided that 
the parameter m is chosen such that Am = Gtr -|- 1. In the presence of boundaries, we derive boundary conditions 
in section (0 that make sure that in this case the resulting IBVP is well posed. In section ljVI|) we show that the 
initial value problem with the gauge conditions (b) is well posed provided that some specified conditions on m and 
the functions /, G and H are satisfied. Symmetric hyperbolic first order formulations of Einstein's equations that 
incorporate gauge conditions that are similar to (b) have been worked out in [29j. 



IV. FIRST ORDER SYMMETRIC HYPERBOLIC FORM (FROZEN SHIFT) 



In this section we recast the BSSN equations with the gauge conditions (a) into a first order symmetric hyperbohc 
system. In order to do so we introduce the extra variables 







dk = I2dk(t>, dkij 


and rewrite Eos. (|5I6I7|I as 


doK = -ae-'^'^f^iAj+l.o., 


doA.j = ae-'^''' 


-i7''5,4,+C7''cf 


dot' = 2a 
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(m l)dkA''' D'K 
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dkji] , Ak 



a 



TF 



+ 1.0.. 



(24) 

(25) 
(26) 

(27) 



where l.o. refers to lower order terms that depend on 0, 7y, K, Aij, F', a, dk, dkij, Ak but not their derivatives. 
Here, we have added the constraint variables Cf^,^- = d[idk]ij with an arbitrary parameter (^ in the equation for Aij. 
As we will see shortly, the addition of these constraints will allow us to obtain a larger family of symmetric hyperbolic 
formulations. Evolution equations for the extra variables are obtained by applying the operator do on the equations 
(Pljl . using the commutation relation [do,dk] — dkP^ ■ di and using the evolution equations H3I4I19|1 for 0, 7ij and a. 
The result is 



dodk = ~2a{dk + Ak)K + didkP' + 2dkdil3\ 

..- _>^2~ 2 

dodkrj = -2a{dk + Ak)A,j + dujdkP' + 2dki{idj)l3^ - -dkt:jdil3^ + 2\,dj)dkl3^ - -7zjdkdiP\ 

doAk - -2aadkK~a—Ak~^-^+AidkP'. 
oa OT 



(28) 
(29) 

(30) 



We have rewritten the BSSN equations (with a fixed prescribed shift but a live condition for the lapse) as a first 
order quasi-linear evolution system for the variables u — {(p,'^ij,a, K, Aij,r'' ,dk,dkij, Ak)^ which is given by the 
equations (|3I4I19I25I26I27I I28I29I30|| . It has the form 



Squ — aA^{u)diU + F{u), 



(31) 



where the matrix-valued functions A*(u), i = 1,2,3, and the vector-valued function F{u) depend on u but not their 
derivatives. An important point to notice here is that we have not added any of the constraints variables H8I10|I to the 
right-hand side (RHS) of the evolution equations for the extra variables. As a consequence, the additional constraints, 
defined by. 



ct 



dk - 12dk^ = 0, 



Ckrj = dktj - dkii-j = 0, 



■'kij 

Ct 



Ak 



a 







(32) 
(33) 
(34) 



that arise when writing the system as a first order one propagate independently on whether or not the remaining 
constraints are satisfied: 



doCi = -2aKC^ + Cfdkf3^ 



doCk. 



kij 

"(A 
■0>^k 



doC^ 



2aA,Ci + Ci^dkH' + 2C,^;(,a,)/3' - -Ct,dil5\ 
OF 



da " 



CfdkP' 



(35) 
(36) 

(37) 



This means that if initial data is given such that C^ — 0, C^^, — 0, C^ — (and suitable boundary conditions are 
chosen), these constraints will also be satisfied at later times and we obtain a solution of the BSSN equations (|3I7|I . 
This is true even if the initial data violates the constraints H = 0, Mi = 0, Cp = of the BSSN system. 



Having obtained a first order quasi-linear system that yields the same solutions than the BSSN system (provided 
that the constraints Cf = 0, C^^ ■ = 0, C^ = are satisfied initially) we now analyze for what range of the parameters 
m, a and C the first order system is symmetric hyperbolic. Introducing the principal symbol A(n) = A^rii where 
n = Ukdx^ is any normalized one-form, this means that we have to find a positive definite matrix H = H(u,a;'') 
which depends smoothly on u and the spacetime coordinates x'^ such that HA(n) is symmetric for all u, cc^ and 
all normalized one-forms nJ33- A necessary condition for this is that each A(n) is diagonalizable and has only real 
eigenvalues. So we first analyze the eigenvalue problem 



fiu — A(n)u. 



(38) 



Explicitly, we have 



II(f> 


= 0, 


IJ-li] 


= 0, 


fia 
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fidk 
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4:171 
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/if, = 2{m^l)e^^An^-—n,K, 
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-ZrikAij , 
-2ankK, 
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'-ntid'^ 



(»"j)fc 
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(39) 
(40) 
(41) 
(42) 

(43) 

(44) 

(45) 
(46) 
(47) 



where An = AiU^, dnij — dkijn^ etc. and f^ = 7yT''. Here, and in the following, we normalize rii with respect to 
the three metric -fij . A convenient way for obtaining the nonzero eigenvalues is by deriving a closed equation for the 
extrinsic curvature. Introducing Kij ~ e'^'^Aij + ^ijK/3 we obtain 



/i Kij = Kij + 2{m — l)nfiKj\„ + (1 — 2m + 2a)ninjK + —{m — l)^ij{K 



K„, 



(48) 



In 24] it was shown that the system is strongly hyperbolic if the operator on the RHS is diagonalizable and has only 
strictly positive eigenvalues. This is the case if and only if the squares of the eigenspeeds. 



fJ'i = 2cr, ^J.2 = 



Am -I 



Ms = "^> 



Ai4 = 1, 



(49) 



are strictly positive, that is, if and only if ttt, > 1/4 and cr > 0. Notice that these conditions are independent of C and 
that for a = 1/2 (that is, if F = aK -\- Fo(a, x^)) and m = 1 all speeds are one or zero. 

In order to find the most general symmetrizer it is convenient to define Kij = e'^'^Aij, to decompose 



3 1 

dkij = -2e^^^ekij + Tlk(ibj) - -Jijbk , 



~ki 



T'dk 



y ' 



where ekij is completely trace-free, and to replace F^, di, bi by the combinations 



= r.j - -di 
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— r, — di 

6 

— adi — Ai 



A,- 



9C + 6, 



20 

A^ + (m - 1)6, , 



Here, we assume that cr > and that 20?n -I- 9<^ — 14 > which implies that the transformation is regular. The first 
condition is necessary for strong hyperbolicity, and the second one can always be achieved by choosing the parameter 
C (which does not appear in the original BSSN system) to be sufficiently large. 
In terms of these variables the non-trivial block of the principal part reads 



liK 



-An 



(50) 



liA^ = -2aniK, (51) 

fikij = Cmj - Ce{ij)n + [n{iVj)] , (52) 

3 - 1 - 

^Ckij = UkK^j - -ji,{^Kj'^„ + -jijKkn , (53) 

/ 9C- 14\ ~ / 1-4to\ , , 

liv^ = (2m+ ^^^ jKm+ (2g+ — \n,K, (54) 

(and ^(f> = 0, ^J'Jij — 0, ^a — 0, liWi — 0, fizi — 0). From this representation of the principal part it is not difficult 
to see that the system is symmetric hyperbolic if and only if 

4TO==6cr+l, cr>0, (55) 

and that in this case a symmetrizer H = H(7'^ , a, m, C) is given by 



(«(i))^Hi.(2) ^ 0(1)0(2) + f S^-'^W^^f + a(i)a(2) + Y^wi'^wf^ + i^' z^hf + 2aK^'^K^'^ + Y^ A^^ Af 



+ y^y'^^^if + 7^-7- (egel^l - Cege^^]) + (2m + 
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7-'u,> "yj 



In order for H to be positive definite we need — 2 < C < 1. (This can be seen by using the orthogonal decomposition 
Ckij = e|jj- +e'^ij, where ef.^j = e(kij) is totally symmetric, and by noticing that e?.-sj, = —el^j/2.) Therefore, we have 
to choose 

max|-2,l-^| <C< 1. (56) 

Since a > this choice is always possible. Summarizing, we have shown that our first order system is symmetric 
hyperbolic if Am — 6cr + 1 > 1, (^ satisfies the inequality (|55|l and if a and C depend smoothly on u and the spacetime 
coordinates x^. This implies that in those cases the corresponding initial value problem is well posed. Since the 
additional constraints propagate, the same result holds for the BSSN system with the gauge conditions (a) when 
Am = 6ct + 1 > 1 and a depends smoothly on u and a;^. Since in this case the evolution system for the constraint 
variables can be reduced to a symmetric hyperbolic system, it follows that the constraints are satisfied if satisfied 
initially. Notice that if < cr < 1/2, there are no superluminal speeds. In the next section, we assume the presence 
of artificial boundaries and discuss boundary conditions. 

V. BOUNDARY CONDITIONS 

Consider the BSSN system (|3I7|) on a bounded domain fi C M'^ with smooth boundary dH.. Consider the slicing 
condition (|19|l with cr > and choose m such that 4m = 6cr + 1. We also assume that the shift is a priori specified, 
and that at the boundary, the shift is tangential to dH.. 

From the previous section we know that the BSSN system can be reduced to a first order symmetric hyperbolic 
system. For such a system, the specification of maximally dissipative 38] boundary conditions yields a well posed 
initial-boundary value formulation |l3|. Maximally dissipative boundary conditions consist in a coupling of the 
ingoing to the outgoing characteristic fields with respect to the normal n to the boundary and some free boundary 
data. The characteristic fields with respect to the normal to the boundary are defined as the projections of u onto 
the corresponding eigenspaces of A(n). 

In order to find the characteristic fields for our first order system, we define 



Ei-i = Cnij - Ce(.y> + ["(jWj)] 



TF 



1 - C r ~ iTF 



"(ij)n 



"C'^j) - gn(idj) - n(,^j) - i^n^^.d^.^^ 
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which is trace-free. Equations (|52l 1531 I54|l imply that 

^ik,j = Eij , (57) 

^iE^j = kij + 2(m - 1) ( n(^,kj)n - ;r7y-^rm ] • (58) 
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In terms of a triad ei, 62, 63 which is such that e\ — n\ it follows from this and Eqs. H50II51|) that the characteristic 
fields with respect to the normal n^ that have nonzero speeds are given by 

l/W = KT^^l^Ar,, (59) 

Vit^ = knn±l^2^E^n. (60) 

K^A' = K^A ± M3"'£^nA , (61) 

Kab ± Eab 



y(±) 



tf 

(62) 



where A,B refer to the triad indices 2 and 3, where [...]*^ denotes the trace-free part with respect to the two- 
dimensional metric 6ab and where /ii, ^2, A*3 are given by the positive square roots of the expressions in (|49|l . A 
short calculation shows that 

u^HA{n)u = ,/2^.((F(+))^-(0-))^)+^ ({V^^^r-ivL^r 



+ M3 s^^ v;(+V,:+) - vtA^vtB') + -^ s^-s-- vd+VW - vjcjvtj) (63) 



1 
2 

The maximally dissipative boundary conditions are given as follows: Let p e 90, and let n^ be the unit outward 
normal to dfl. Then, the boundary conditions at p are 

(64) 
(65) 
(66) 
(67) 

where a, b are smaller than one in magnitude and the matrices c^ and (Iab have norm smaller than one, and where 
G, Gnn, GnA and Gab are freely specified source functions (subject to the condition S^^Gab — 0). In order to 
illustrate why these boundary conditions lead to a well posed IBVP, let us linearize the equations around an arbitrary 
background. The resulting equations have the form 

Oqv = aA^diV 4- Bw, 
where v denotes the perturbation. Defining the energy norm 



T/(+) = 


= ay(")-|-G, 


^ nn 


= bV,[~^ + G„„ , 


^nA - 


= '^IVnB + GnA , 


y{+) _ 

^ab - 


- d'ily^cji^GAB 



In 
taking a time derivative, using the symmetries of the matrices H and HA* and using Gauss' theorem, we find 

^£ ^ 2 [ v'^UHaA' + (3')diV + Bv] d^x 
dt Jn 

= I [d, {v^ailA'v + v^np'v) + v^ (HB + B^H - a,(aHA^ + H/^')) v] (fix 
Jn 

< [ av'^TiA{n)vd^x + CE, (68) 

Jdn 

where we have used the fact that the shift is tangential to the boundary at the boundary and where G is a constant 
that only depends on bounds for B and H^^5i(aHA' -I- H/3*). If the boundary conditions are homogeneous, i.e. if 
G = 0, Gnn = 0, GnA = 0, Gab = 0, Eqs. (|63I64I65I66I67|) immediately imply that the boundary integral is negative 
or zero, and we obtain the energy estimate £{t) < exp(Gi)£(0). If the boundary conditions are inhomogeneous one 
can bound £{t) by £{0) and the L2-norm of the boundary data [13, y^l- These energy estimates play a key role 
in proofs for well posedness. These proofs can be generalized to quasi-linear symmetric hyperbolic systems, see for 
instance p^ . 

Therefore, the boundary conditions (|64I65I66I67() lead to a well posed initial-boundary value formulation. Since 
the shift is tangential to the boundary at dfl, the additional constraints propagate as before, and thus the same 
boundary conditions applied to the equations (|3I7|I . where we perform the replacements (|24ll . yields a well posed 
initial-boundary formulation for the BSSN system. In particular, choosing a = 6 = 0, c^ = 0, dAB — 0? and setting 



the source functions G, G„„, GnA, Gab to zero, corresponds to Sommerfeld-type boundary conditions, in the sense 
that these conditions are algebraic conditions for the first order systems which are perfectly absorbing for plane waves 
of normal incidence to the boundary in the frozen coefficient approximation. Explicitly, we obtain the six boundary 
conditions 

K ^- n'd^a = 0, 



5'^^ A: 



2a a 
V4m-1 "' 
'm 



e^eij - ^ bAB^'^^'eh^'n 



ij + ^ij 



= 0, 
= 0, 
= 0, 



A = 2,3, 
A,B = 2,3, 



(69) 

(70) 

(71) 
(72) 



where 



E,j^--e^'^n''dk^,,+ 



1 c 



TF 



where n' is the unit outward normal to the boundary and the vectors 62 and 63 must be chosen such that n^, e\, e\ 
form a triad with respect to the three metric 7.^, and rii = "fijU^ . The vectors 62 and 63 are unique up to a rotation; 
such a rotation does not alter the boundary conditions. The parameter ( has to be chosen such that the inequality 
H56|l is satisfied. The boundary conditions (|69I7()I71I72(I can be generalized to inhomogeneous conditions by replacing 
the zeroes on their right-hand sides by freely specifiable source functions G, G„„, G„yi, Gab- If the solution is known 
in a neighborhood of the boundary, one can compute these source functions by evaluating the left-hand sides of Eqs. 
(|69I70I71I72() . Notice that the occurrence of the parameter (, which does not appear in the BSSN system, in the 
boundary conditions has its origin in the ^-dependence of the unphysical energy £ defined by the symmetrizer. 

VI. STRONG HYPERBOLICITY WITH A DYNAMICAL SHIFT 

Here we consider the BSSN equations (|3I7|I with the live gauge conditions (b) , see section IIIII In this case one 
could proceed as in the frozen shift case and introduce the shift and its first derivatives (with respect to time and 
space) as extra variables. One obtains a first order system that is equivalent to the original system provided that the 
additional constraints arc satisfied. Unfortunately, we did not succeed in finding a symmetrizer for the resulting first 
order system. Our goal in this section, therefore, is more modest: We show that the BSSN system with the live gauge 
conditions is strongly hyperbolic and so prove that the resulting Cauchy problem (in the absence of boundaries) is 
well posed. A related analysis for a different form of the system has been performed in [33| . 

For differential equations that are not first order, a definition of strong hyperbolicity has recently be given in 19, 22] 
that does not require the introduction of extra variables (nor extra constraints). It is based on pseudo-differential 
calculus. The intuitive idea behind this definition is to freeze the coefficients in the differential equations at some 
fixed point and to analyze the resulting linear constant coefficient problem by means of a Fourier transformation in 
space. In our case, the frozen coefficient problem is given by 



do(i> 
doK 

doArs 



-6^+6-^^ 



2i 



-2aArs + 2i^k{r^s)P'' - —Irst^kP 



- ^.-40 



-,fei, 



-7" UJkUJars +ilk(r^s)^' + 2u}rU}s4> + ^r^ s — 



TF 



+ 1.0., 



doa — 



-^'^^ujkLUiP" - -Y'uj^ujkP'' + 2a 



2im , 



i{m - l)ujkA'"' Y'^rK 



+ 1.0.. 



a^G{a,(f>,x'')B', 



2im' 

i{m' - l)uJkA^-' - —-Y'iOrK 
o 



+ 1.0., 
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where a hat denotes the Fourier transformation in space, (j){io) = J (t>{x) exp(— iw • x) (Px, and l.o. denotes terms that 
depend on lower order spatial derivatives. Here, we have also allowed for a parameter m' that is different than m 
in the evolution equation for B^. We can rewrite this as a first order system in t and LUi by writing LUi = It^jrii, 
\lu\ — yJ^^^bJk'^i, and introducing the variables 



V 



/I — 



bs - 

f. = 



i\uj\ 



^ Irs 5 



la \uj\a, 



Bs - 7r.S^ 

In terms of the variables we obtain a first order pseudo-differential system of the form 

dtu — i\uj\ (aP(n) + j3^ni) u + l.o., 

where ii = {(f, hij, K, kij, a, hi, Ti, Bif" . The system is strongly hyperbolic if there exists a positive definite Hermitian 
matrix H(a;^, u, n) which is smooth in all its entries such that HP is symmetric. A necessary condition for this is that 
P is diagonalizable and has only real eigenvalues. Therefore, we first consider the eigenvalue problem /iu — P(n)u; 
explicitly 



liip 



■K 



fX firs — l^rs 

II K = —a, 

jX f^rs f^rs 

fia = —fK, 

libs = GBs, 



l{rOs) 



TF 



'^(r-Ts) ^ 2nrns(p — nrUgil 



TF 



1 



bs + -r^s^n + 2(m - l)k„ 



Am 



-n^K, 



liBs = H 



1 - - 4r7j' 
bs + ■:;nsbn + 2(m' - l)kns :T-nsK 



where 6„ — j^'^rirbs and fc„j — ^^'^rirksj. A careful analysis reveals that the matrix on the RHS has the eigenvalues 
0, ±^1, ±/i2, ±y(i3, ±^4, ±/i5 where 

Ml = v//, 



M2 



Am — 1 



M3 



/i4 



M5 



'GH 



M6 



AGH 



Therefore, we need m > 1/4, / > and GH > 0. {Ii G = H = the equation for the shift decouples, and we are 
back in the case considered in the previous section). Furthermore, it turns out that the matrix is diagonalizable only 
if AGH y^ 3/ and provided that m' = 1 if vn = GH or AGH = Am — 1. In the remaining cases the system is only 
weakly hyperbolic which, in the nonlinear case, can lead to exponential growth with arbitrarily small growth time. 
Introducing the factors 



^3 - 



AGH 



3/ - AGH ' 

Gjm' - 1) 
4m - 1 - AGH ' 
2(m' - 1)GH 

m-GH ' 



if Am — 1 7^ AGH and Q2 arbitrary otherwise, 
if TO 7^ GH and CI3 arbitrary otherwise, 
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the eigenfields can be expressed as 

Zo = Sniip — 2{m — l)hr, 



Z, = H 



V^' 



V. 



(±) 



2{m — m')h„i + m'Ti 



mBi , 



KTf^i^a , 



y{±) _ L ^^^..-iff, _2^ , 1^ 
^nA = ^nA T l^-^^ ( Ka - ^f A 



y(±) _ 



kAB T M4 'iAS 



*/ 



bA - ^zknA ± /^5 ^ 



GBa + r!3 ( Ka - ^Va 



2K 



yi±) = fe„ + riiX - r!2 4m - — ± M^' 



2 ~ 4 

GS„ - Vtia + ri2 ( hnn - ^r„ + -1^ 



where the components n, A — 2, 3, refer to triad indices as before. The matrix H which symmetrizes P can now be 
built by summing over the square of the eigenfields. It is smooth in the one- form n^. In order to see this one uses, for 
example, hnAh^^ = n*n'^(7-'' — n^n^)hijhki. It is also smooth in the other variables provided that f^i, ^2 and fl^ can 
be chosen smoothly and are bounded. A simple possibility of achieving this is by choosing m = ni' — I and / — kGH 
with K a constant that is unequal 4/3. The pseudo-differential calculus shows that in these cases the full nonlinear 
Cauchy problem is well posed. Since the evolution system for the constraint variables can be reduced to a symmetric 
hyperbolic system ii m > 1/4 it follows that the constraints are satisfied if satisfied initially. 



VII. CONCLUSION 



We discussed some mathematical aspects of the BSSN system which is currently used by several groups in numerical 
relativity. In particular, we derived a well posed initial-boundary value formulation of the BSSN system with a 
Bona-Masso like slicing condition for the lapse and a frozen shift. This is achieved by introducing extra variables and 
recasting the evolution equations into a first order symmetric hyperbolic system, for which maximally dissipative 
boundary conditions are specified. The introduction of extra fields brings additional constraints, and the original 
BSSN system and the first order symmetric hyperbolic system derived in this article are only equivalent if these 
constraints are satisfied. However, we showed that the associated constraint variables obey a closed evolution system 
that is independent of the other constraints. Moreover, by choosing the shift to be tangential to the boundary, these 
additional constraints propagate tangentially to the boundary. This implies that they are satisfied everywhere at 
later times if satisfied initially, even if the other constraints are violated. This allows us to return to the second order 
system and to conclude that the BSSN system with the specified boundary conditions is well posed; in particular, 
unique solutions local in time exist, and depend continuously on the initial and boundary data. To our knowledge, 
the specified (six) boundary conditions H69I70I71I72() have not yet appeared in the literature. 

In general, the boundary conditions derived in this article are not compatible with the constraints of the 
BSSN system. They can feed in some constraint violating modes. Nevertheless, they are consistent with the evolution 
equations and constitute a first step towards improving numerical evolutions of the BSSN system. In particular, the 
present analysis offers the possibilityto construct constraint-preserving boundary conditions [Ijl in the linearized 
case, following the lines of \iaL Il6l Il7l |. Furthermore, the derivation of the symmetrizer and the energy estimate 



presented in section |V| should be useful as a guidance 
numerical stability at least at the linearized level 



rinciple to construct discretizations schemes that guarantee 



ance principle 

liillillliil 



We have also considered dynamical gauge conditions for lapse and shift and obtained a class of second order 
evolution equations which can be shown to be strongly hyperbolic using pseudo-differential calculus. For these 
systems, one can show well posedness of the initial value problem. Here, the presence of boundaries has not been 
considered. To derive boundary conditions in the case of a finite domain one could proceed as follows: First, derive 
a first order system by introducing extra variables as described at the beginning of section IVll Next, consider the 
matrix A(n) multiplying the derivatives normal to the boundaries. In case this matrix is diagonalizable, to every 
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strictly positive eigenvalue of A(n) there corresponds a Sommerfeld-type outgoing boundary condition given by the 
condition of a vanishing projection of the field vector onto the corresponding eigenspace. This corresponds to setting 
zero the incoming characteristic fields with respect to the direction which is normal to the boundary. Finally, well 
posedness of the initial-boundary value problem in a suitable Hilbert space has to be proved J39|. Necessary conditions 
for well posedness can be obtained by using the method of Laplace transformation, see, for example, ;10i. 30J . The 
derivation of boundary conditions in the dynamical shift case is beyond the scope of the present work. 

The gauge conditions considered here differ from the ones used in j21j for numerical simulations only by the 
replacement dt f—^ dt — (3^ dj , which leads to a simpler principal part and makes it more amendable to analyze the 
algebraic conditions that guarantee symmetric or strong hyperbolicity. Preliminary investigations of the problem 
without this replacement have been done in |35j where one of the Sommerfeld-type conditions has already been 
computed. The structure of this condition is more complicated than the conditions derived in this article. 
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